Velocity of front propagation in the epidemic model A + B — > 2A. 
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We study front propagation in the irreversible epidemic model A + B — > 2A in one dimension . 
Here, we allow the particles A and B to diffuse with rates Da and Db, which, in general, may be 
different. We find analytic estimates for the front velocity by writing truncated master equation in a 
frame moving with the rightmost A particle. The results obtained are in reasonable agreement with 
the simulation results and are amenable to systematic improvement. We also observe a crossover 
from the linear dependence of front velocity V on Da for smaller values of Da to V oc \/Da for 
larger Da, but numerically still significantly different from the mean field value. The deviations 
reflect the role of internal fluctuations which is neglected in the mean field description. 



I. INTRODUCTION 

In many natural phenomena we encounter propagat- 
ing fronts separating different phases 0. In this paper 
we study the autocatalytic reaction A + B — ► 2A in a 
one dimensional lattice. Models of this type have been 
studied activelybecause of several interesting features of 
the process 0, 0, ■ We can think of this reaction as 
a model for various spreading phenomena such as an epi- 
demic: particle A is an infectious agent which infects a 
healthy particle B on contact. We are interested in find- 
ing the velocity of propgagation of this infection. The 
macroscopic description for this process is given by the 
following mean field equation 0. 



dpA 
dt 



D A Ap A + ep A p B = D A Ap A + ep A (p - p A ) (1) 



Here, p A (x,t) is the local density of A particle at posi- 
tion x and time t, D A is its diffusion coefficient and e is 
the reaction rate of the particle. In the second line we 
have used p A + pb = Po which is constant. Equation(2J 
arises in the macroscopic description of diverse physical 
processes and serves as a generic model for front prop- 
agation in a system in which a stable state invades an 
unstable state 0. Equation (JJJ allows a family of trav- 
elling wave solutions with speed V > 2^/ep^D A invading 
the unstable state p A = from the stable state p A = po- 
However, for steep enough initial conditions, the selected 
asymptotic speed is the minimum speed Vq = 2y / epoD A 
0. In the microscopic lattice model in low dimensions, 
the discreteness effects alter the dynamical properties of 
the propagating front and the results obtained are dif- 
ferent from that using deterministic mean field equation 
00. Here, in this paper we will present the Monte 
Carlo simulation results and approximate analytic pre- 
diction for the front velocity in a one dimensional lattice 
model for A + B — > 2 A defined in the next section (II). 



II. MODEL AND FRONT VELOCITY 

In our simulation, we start with a ID lattice [0, oo) ini- 
tially randomly filled with B particles with overall con- 



centration pq. We have taken hard core exclusion into 
account: any site can be occupied by maximum one par- 
ticle. At the left end of the lattice we place an A par- 
ticle. We allow both the particles A and B to diffuse 
to the nearest neighbour empty sites with rates D A and 
Db respectively. When a 'sick' particle A encounters the 
'healthy' one B sitting on a neighbouring site, then B 
gets infected with rate e and is converted to an A parti- 
cle. In this way the front ( rightmost A ) propagates and 
we wish to find the velocity of this front. 

In Ref 0, authors studied a similar problem using 
Smouluchowski approach and obtained the steady 
state density for B and A particles at a distance x from 
the moving front particle (at x — 0) as: pb(x > 0) — 
p (l - e - Vx /( D A+D B )} anc j pA ( x < o) — pA respectively. 
By invoking the equality of particle fluxes to the left (j _ = 
V PA ) and right[(j+ = (D A + D B )(d/dx)p B (+0))} they 
found p A — po, but it does not fix the front propagation 
velocity V. In the present paper, we try to find the ex- 
pression for the front velocity by visualising the front as 
random walker moving with forward and backward rates 
as P + and P_ respectively, and thus speed of the front 
V = P+ — P- . In order to compute P+ , we notice that 
there are two ways in which the front will move to the 
right: (1) when the front diffuses to the right neighbour- 
ing site, provided it is empty with rate D A (1 — p^), here 
pi is the probability that the site just ahead of the front is 
occupied, (2) when the site just ahead of front is occupied 
by B particle, then due to reaction it changes the neigh- 
bouring B into A and thus the front moves forward with 
rate ep\ and hence P+ = ep\ + D A (1 — pf). Similarly, 
the front moves backward due to its diffusive move to the 
left empty neighbouring site and hence P_ = D A (1 — p^), 
p\ is the probability that site just behind the front is oc- 
cupied. Using these values for P + and P_, we have for 
the asymptotic velocity, 



V = ep+ + D A (p^ -p+), 



(2) 



which is exact, given pf and p^. The multi-particle na- 
ture of the problem is embodied in finding p\ and p^ . 
The fact that equation J2J) is exact can be seen from 
the excellent agreement between simulation result of V 
and its estimate using values of p\ and p^ obtained di- 
rectly from simulation, Fig.Q). Here, we also observe a 



2 



0.16 
0.12 
> 0.08 
0.04 

t 





r '0.02 
0.01 



« » . . 



0.005 0.155 0.305 



D 



A 



0.1 



0.2 



0.3 



0.4 



FIG. 1: Front velocity V as a function of Da for Db = 
0.05, e = 0.20, po = 0.5. The broken line represents the front 
velocity obtained by simulation while the filled circles corre- 
sponds to the velocity obtained by using Eq. J^J and substi- 
tuting the values of pf and from simulation. Inset: The 
differenece \p± — po as a function of Da, p\ taken directly 
from simulation. 
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FIG. 2: Log-log plot for the front velocity as a function of 
D A for D B = 0.05, e = 0.20, p = 0.5. Here V ~ Da and we 
notice a crossover from a ~ 1.0 to a m 0.50. 



crossover from linear dependence of V on Da for small 
Da to V oc for large as shown in Fig. 

although we have verified that V 7^ Vb(mean field value). 

From 0, the density profile of A particles behind the 
front is constant and equal to po- Thus we can make 
approximation p\ = po. But we note that this approxi- 
mation fails to work when Da = where we get pa = 1 at 
all the sites behind the front irrespective of po. Thus the 
approximation p^ = po works better for larger values of 
Da while for smaller Da it appears to be a poor approx- 
imation as shown in the Fig. Q inset). Now using the 
approximation p^ = po we try to find analytic expression 
for p\ . In the following we find an approximate ana lyti c 
expression for p\ using the technique developed in [ldj . 
The technique involves writing a truncated master equa- 
tion in a frame moving with the front ( i. e. the rightmost 
A particle). For example, the simplest one is the set of 
two states {Acj>,AB}. In this truncated representation 
each state contains the occupancy at two sites (I = 2), 
the leftmost A represents the front particle while <p and 
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FIG. 3: Microscopic moves: Filled circles/squares represent 
A/B particles, (a) Diffusion of A to the left empty site lead- 
ing to transition from the state {AB) to (A<f>), (b) Diffusive 
jump of B particle to its right empty site leading to transition 
from (AB) to (A<j>) , (c) Transition from (AB) to (A<j>) due to 
reaction and (d) Diffusion of B to its left and thus making 
transition from (A(j>) to (AB). 



B stand for empty and a site occupied by B just ahead 
of the front respectively. These two states make transi- 
tions between each other due to reaction or diffusion of 
particles, see Fig. ©. For example, in the transition 
shown in Fig. (|3Ji), the configuration AB changes to A<j> 
if A particle diffuses to its left empty site and it takes 
place with rate Da(1 — Pi)- In Fig- El 3 ), diffusion of B 
particle in the configuration AB to its right neighbouring 
empty site changes AB to A<f>. If the probability of oc- 
cupancy of B particle at a second site ahead of the front 
is denoted by p\, then AB — > A(f> with rate £>s(l — p\). 
In the Fig. |(Ht)) when A infects B in the configuration 
AB then AB — > A<j) if the second site ahead of the front 
is empty and this occurs with rate e(l — p^)- Similarly, 
in Fig(|3|i), Atfi changes to AB with rate Dbp\- Thus 
considering all the transitions between the two states A<p 
and AB and approximating — po, the following mas- 
ter equation may be written [ll|: 



= {D B (1 - p+) + D A (1 - po) + e(l - pi)}p AB 
- (D A + D B )ptPA<f>, 

= (D A + D B )p+ P A4, - {D B (l - p+) + D A (l - Po) 
+ e(l-p+)} PAB . (3) 



df 



dpAB 

dt 



In the steady state, by definition, p\ = pab and thus 
from Eq. @ with the normalization pa$ + Pab = 1 we 
obtain 



Pi = 



(D A + D B )p\ 



{D B + D A pt + D A (l - po) + e(l - pj)} 



(4) 



Thus, in order to find p\ we need to know and as 
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a first approximation assuming p\ = po, we find 
(D A + D B ) Po 



Pi 



D A + D B + e(l - po) 



(5) 



Using this approximation for pf in Eq. we obtain the 
following expression for the velocity, 



V = 



ep (2D A + D B -D A p ) 
D A + D B + e(l - po) 



(6) 



We have shown the estimate for the velocity as obtained 
from Eq. in Fig. J1J. Here we note that when e = D A , 
the above approximate result for V satisfies the simula- 
tion result nicely. However, as the magnitude of e — D A 
increases the analytic prediction for V departs gradually 
from the simulation results. From Eq. we note that 
when D A — e the expression for V is independent of pj 
and hence we obatain nice agreement between the theory 
and simulation results. But the dependece of V on pf in- 
creases gradually with the inceasing magnitute of e — D A 
and hence we observe poor agreement between analytic 
and simulation results as we move away from the point 
D A = e. In order to get better prediction for the ve- 
locity, especially for larger values of e — D A , we need to 
study the states having larger number of particles. For 
example, we studied the following set of states(7 = 3): 
{A4>4>,AB4>,A4>B,ABB}. If we denote the probability 
of occupancy at third site ahead of the front as pi , the 
evolution of these states is described by the following 
master equation. 



{D B (l - p+) + D A (1 - p )}p M B 
e (! - pt)PAB<t> - (D A + D B )p^p A4>4> , 
(D A + D B )p+p Mlt> + {D B + D A (1 - po) 
ep^}p ABl p + D A (1 - pq)pabb - {2D A 



dt 
dp Al p B 



dt 



+ 2D B - D A p - D B pj}p A<j)B , 

dj = {D B + D A {1- p+)}p MB + {D B {l- p+) 

+ e(l - pj)}p ABB - {D A + D B + e 

- D A p + D B p% }p AB<j> , 



dp ABB 
dt 



D A p£p Alj>B + D B p^p AB4> - {D A (l - p ) 



D B (1- p+) + e(l- p+)}p ABB . 



(7) 



In order to solve the Eqs. we need to know pi . To 
find the steady state solution of Eq. 10, and thus pi = 
PAB4> +Pabb, one could approximate pi = po, as was 
done for I = 2 case. However, a better estimate for p^ 
can be obtained by writing a mean field equation for the 
densities an d P3 at second and third site ahead of the 



front, 



d 4r = ( P i+pl)(l-pl)D B -D B p+(2-p+-p+) 7 



dt 
dt 



do + 

= (pi + pi)(i- pi)D B - d bP +(2- p+ - p+). 

(8) 



In steady state Eqs. JSJ) gives, p^ = {pi + 2p±)/3. Now 
using the value of pf from Eq. (J3J) and approximating 
Pa = Po, w e obtain the following approximate expression 
for p£. 

+ _ SDaPo + 3D B p a + 2ep (l - p Q ) 
P3 ~ 3{D A + D B + e(l- PQ )} U 

Using this value of pi in Eqs. (J2J), we obtain the steady 
state solution satisfying the normalisation pA4>4> +PA41B + 
PAB<f> +Pabb = 1< Once wc know these probabilities the 
density at site just ahead of the front is obtained as: 
pi = PABip +Pabb- Substituting this value of pi in Eq. 
P)l. we obtained an improved analytic estimates for the 
velocity as shown in the Fig. Q. Here we notice that as I 
increases the data corresponding to the analytic estimate 
comes closer to the simulation result. In this approach 
one can systematically improve upon the estimates for 
pi and thus V by including states having successively 
larger number of particles. We have also observed that 
when Da=0, the front position x grows with time t as 
t x l 2 as was reported in Q. Hence, when Da = 0, the 
asymptotic velocity V approaches to zero. Using Eq. (J2J 
for D A = 0, velocity is given as V — epi . The condition 
for zero velocity simply implies that pi — > 0. In the Fig. 
|(5J| , we have shown the simulation results for the density 
profile ahead of the front for different values of D A and 
notice that as Da — * 0, pi — > 0, which gives V — > 0. We 
also note that when Da = 0, A particles form a compact 
cluster while B particles keep on moving diffusively. The 
number of A particles created during time interval dt, 
dNA = Vdt should be equal to the number of consumed 
B particles, dN B = Vpodt. However, the equality dN A = 
dN B holds for non zero V only if po = 1 while for po 7^ 1 , 
it gives V = 0. We also notice that when Da = 0, the 
front dynamics is like a DLA in Id, where, A particles 
are static and growth occurs when B particle diffuses and 
sticks with the neighbouring A particle. 



III. CONCLUSION 

In this paper we have analyzed the reaction A + B — > 
2A using an approximate analytic method and Monte 
Carlo simulation. Considering the front as the position 
of rightmost A particle we have been able to obtain ana- 
lytic estimate for its velocity by writing master equation 
in the frame moving with the front. One can improve 
upon the results systematically by including states with 
successively larger number of sites in front of the leading 
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A particle. The mean field equation gives V oc y/D~A, 
but we observe linear dependence, V oc Da for small 
Da reflecting the important role of discrete fluctuations 
which is neglected in the deterministic mean field equa- 
tion 0. 



0.1 0.2 0.3 0.4 

Da 

FIG. 4: Velocity V as a function of Da for Db = 0.05, e = 
0.20, po = 0.5. Filled circle(») represents analytic estimate 
for I — 2, that is, for states {A(j>, AB} . Filled square corre- 
sponds to the results obtained for I = 3, that is, for the states: 
{A(/>(j>, A(j>B, AB<j>, ABB} while the open circle (o) is the sim- 
ulation results for the front velocity. Here, we note that as 
I increases the analytic prediction curve comes closer to the 
simulation results. 
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FIG. 5: Density profile ahead of the front particle as a func- 
tion of Da while keeping Db = 0.05, e = 0.20, and po = 0.5 
fixed. The simulation data from top to bottom corresponds 
to D A = 0.025, 0.015, 0.005, 0.004, 0.002, 0.0005, 0.0. Here we 
note that as Da — ► 0, p\ — > 0. 
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